#pragma once

#include "Coordinate.hpp"
#include "../Resource/Texture/TextureCube.hpp"
#include "../Resource/Texture/Texture2D.hpp"
#include "SH/SHCoeff.hpp"
//#include "SH/SHCoeff4f_SSE.hpp"
#include <Math/Array.hpp>

namespace zzz{
  template <typename T>class Cubemap;

  struct CubemapPos
  {
    CartesianDirCoordf direction;
    float solidangle;
  };
  typedef vector<CubemapPos> CubemapPosList;

  class Cubemapbase
  {
  public:
    static const CubemapPosList& GetPosList(zuint n);
    static const Cubemap<float>& GetSHBaseCubemap(int l,int m,int size);
    static const Cubemap<float>* GetSHBaseCubemapPointer(int l,int m,int size);
  private:
    static map<zuint,CubemapPosList*> PosListMap;
  private:
    //SH Projection Helper
    typedef Vector<3,int> SHBaseCubemapDesc;
    static map<SHBaseCubemapDesc,Cubemap<float>*> SHBaseCubemap;
    static float SHBase(int l,int m,float theta,float phi);
    static float K(int l,int m);
    static float P(int l,int m,float x);
  };

  template <typename T>
  class Cubemap : public Cubemapbase
  {
  public:
    Cubemap()
    {
      cubesize=0;
      facesize=0;
      datasize=0;
      v=NULL;
    }
    Cubemap(int n)
    {
      cubesize=n;
      facesize=n*n;
      datasize=facesize*6;
      v=new T[cubesize*cubesize*6];
    }
    Cubemap(const Cubemap<T> &cube)
    {
      v=0;
      *this=cube;
    }
    ~Cubemap()
    {
      Clear();
    }
    void Dump()
    {
      for (zuint i=0; i<datasize; i++)
      {
        if (i%facesize==0) zout<<"-------------------------------------------------\n";
        if (i%cubesize==0) zout<<"| ";
        zout<<v[i]<<' ';
        if (i%cubesize==cubesize-1) zout<<"|\n";
      }
    }
    const Cubemap& operator=(const Cubemap<T> &cube)
    {
      ChangeSize(cube.cubesize);
      memcpy(v,cube.v,sizeof(T)*datasize);
      return *this;
    }
    void ChangeSize(int n)
    {
      if (n==cubesize) return;
      Clear();
      cubesize=n;
      facesize=n*n;
      datasize=facesize*6;
      v=new T[datasize];
    }
    void Clear()
    {
      delete[] v;
      v=0;
      cubesize=0;
      facesize=0;
      datasize=0;
    }
    void Zero(const T x=0)
    {
      for (zuint i=0; i<datasize; i++) v[i]=x;
    }
    inline T& At(const CubeCoord &c)
    {
      return v[c.ToIndex()];
    }
    template <typename T1>
    inline T Get(const CubeCoordNormalized<T1> &c) const
    {
      float ratioX,ratioY;
      CubeCoord basePos=c.GetBaseCubeCoord(cubesize,ratioX,ratioY);
      T ret=0;
      CubeCoord pos=basePos;
      pos.Standardize();
      ret+=v[pos.ToIndex()]*(1-ratioX)*(1-ratioY);
      basePos.x++;
      pos=basePos;
      pos.Standardize();
      ret+=v[pos.ToIndex()]*1*(1-ratioY);
      basePos.x--;basePos.y++;
      pos=basePos;
      pos.Standardize();
      ret+=v[pos.ToIndex()]*(1-ratioX)*ratioY;
      basePos.x++;
      pos=basePos;
      pos.Standardize();
      ret+=v[pos.ToIndex()]*ratioX*ratioY;
      return ret;
    }
    inline void Set(const CubeCoordNormalized<T> &c,const T value,bool plus=false)
    {
      float ratioX,ratioY;
      CubeCoord basePos=c.GetBaseCubeCoord(cubesize,ratioX,ratioY);
      T ret=0;
      basePos.Standardize();
      if (plus) v[basePos.ToIndex()]+=value/(4*(1-ratioX)*(1-ratioY));
      else v[basePos.ToIndex()]=value/(4*(1-ratioX)*(1-ratioY));
      basePos.x++;
      basePos.Standardize();
      if (plus) v[basePos.ToIndex()]+=value/(4*1*(1-ratioY));
      else v[basePos.ToIndex()]=value/(4*1*(1-ratioY));
      basePos.x--;basePos.y++;
      basePos.Standardize();
      if (plus) v[basePos.ToIndex()]+=value/(4*(1-ratioX)*ratioY);
      else v[basePos.ToIndex()]=value/(4*(1-ratioX)*ratioY);
      basePos.x++;
      basePos.Standardize();
      if (plus) v[basePos.ToIndex()]+=value/(4*ratioX*ratioY);
      else v[basePos.ToIndex()]=value/(4*ratioX*ratioY);
      return ret;
    }
    inline T& operator[](zuint n)
    {
      return v[n];
    }
    inline const T& operator[](zuint n) const
    {
      return v[n];
    }
    template <typename T1>
    float Dot(const Cubemap<T1> &c) const
    {
      T ret=0;
      for (int i=0; i<datasize; i++)
        ret+=v[i]*c[i];
      return ret;
    }
    bool LoadFromPfm(const string &filename)
    {
      printf("LoadFromPfm: Must be specialized to use\n");
      return false;
    }
    bool SaveToPfm(const string &filename)
    {
      printf("SaveToPfm: Must be specialized to use\n");
      return false;
    }
    void LoadFromTextureCube(TextureCube &t, int channel=0)
    {
      printf("LoadFromTextureCube: Must be specialized to use\n");
    }
    void LoadFromTexture2Ds(Texture2D *t, int channel=0)
    {
      printf("LoadFromTexture2Ds: Must be specialized to use\n");
    }
    void ProjectToHaar()
    {
      Haar2D(v,cubesize);
      Haar2D(v+1*facesize,cubesize);
      Haar2D(v+2*facesize,cubesize);
      Haar2D(v+3*facesize,cubesize);
      Haar2D(v+4*facesize,cubesize);
      Haar2D(v+5*facesize,cubesize);
    }
    void LoadFromHaar()
    {
      AntiHaar2D(v,cubesize);
      AntiHaar2D(v+1*facesize,cubesize);
      AntiHaar2D(v+2*facesize,cubesize);
      AntiHaar2D(v+3*facesize,cubesize);
      AntiHaar2D(v+4*facesize,cubesize);
      AntiHaar2D(v+5*facesize,cubesize);
    }
    void SingleLoadFromHaar(zuint dir)
    {
      memset(v,0,sizeof(T)*datasize);
      zuint face=dir/(facesize);
      dir=dir%(facesize);
      memcpy(v+facesize*face,SingleLoadFromHaarFace(dir),sizeof(T)*facesize);
    }
    T *SingleLoadFromHaarFace(zuint dir)
    {
      static T **pre=NULL;
      static int precubesize=0;
      if (pre==NULL || precubesize!=cubesize)
      {
        zout<<"SingleLoadFromHaarPrecompute!\n";
        if (pre!=NULL)
          for(zuint i=0; i<facesize; i++)
            delete[] pre[i];
        delete[] pre;
        precubesize=cubesize;
        pre=new T*[facesize];
        for (zuint i=0; i<facesize; i++)
        {
          pre[i]=new T[facesize];
          memset(pre[i],0,sizeof(T)*facesize);
          pre[i][i]=1;
          AntiHaar2D(pre[i],cubesize);
        }
      }
      return pre[dir];
    }
    template <int N,typename T1>
    void SingleProjectToSH(int pos, T value, SHCoeff<N,T1> &sh)
    {
      int cur=0;
      for (int l=0; l<N; l++) for (int m=-l; m<=l; m++)
      {
        const Cubemap<float> &shbase=Cubemapbase::GetSHBaseCubemap(l,m,cubesize);
        sh[cur]=value*shbase[pos];
        cur++;
      }
    }
    template <int N,typename T1>
    void ProjectToSH(SHCoeff<N,T1> &sh)
    {
      int cur=0;
      for (int l=0; l<N; l++) for (int m=-l; m<=l; m++)
      {
        const Cubemap<float> &shbase=GetSHBaseCubemap(l,m,cubesize);
        T1 res=0;
        for (zuint i=0; i<datasize; i++)
          res+=v[i]*shbase[i]; 
        sh[cur]=res;
        cur++;
      }
      return;
    }
/*
    void ProjectToSH(SHCoeff4f_SSE &sh)
    {
      int cur=0;
      for (int l=0; l<4; l++) for (int m=-l; m<=l; m++)
      {
        const Cubemapf &shbase=GetSHBaseCubemap(l,m,cubesize);
        float res=0;
        for (zuint i=0; i<datasize; i++)
          res+=v[i]*shbase[i]; 
        sh[cur]=res;
        cur++;
      }
      return;
    }
*/
    template <int N,typename T1>
    void LoadFromSH(SHCoeff<N,T1> &sh)
    {
      Zero();
      int cur=0;
      for (int l=0; l<N; l++) for (int m=-l; m<=l; m++)
      {
        if (sh[cur]!=T1(0))
        {
          const Cubemap<float> &shbase=GetSHBaseCubemap(l,m,cubesize);
          for (zuint i=0; i<datasize; i++)
            v[i]+=sh[cur]*shbase[i]; 
        }
        cur++;
      }
    }
    T *v;
    zuint cubesize, facesize, datasize;
  private:
    //Haar Wavelet Projection Helper
    static void Haar1D(T *vec,int n)
    {
      int w=n;
      T *vecp = new T[n];
      memset(vecp,0,sizeof(int)*n);
      const float sqrt2=sqrt(2.0);
      while(w>1)
      {
        w/=2;
        for(int i=0; i<w; i++)
        {
          vecp[i] = (vec[2*i] + vec[2*i+1])/sqrt2;
          vecp[i+w] = (vec[2*i] - vec[2*i+1])/sqrt2;
        }
        memcpy(vec,vecp,sizeof(T)*w*2);
      }
      delete[] vecp;
    }
    static void Haar2D(T *vec,int n)
    {//2d haar on a nxn map
      int height=n;
      while(height>1)
      {
        for (int i=0; i<height; i++)
          Haar2Dhelper1(vec+n*i,height);
        for (int i=0; i<height; i++)
          Haar2Dhelper2(vec+i,height,n);
        height/=2;
      }
    }
    static void Haar2Dhelper1(T *vec,int n)
    {//one step horizontal haar
      T *vecp=new T[n];
      const float sqrt2=sqrt(2.0);
      const int w=n/2;
      int x=0;
      for (int i=0; i<w; i++)
      {
        vecp[i]=(vec[x]+vec[x+1])/sqrt2;
        vecp[i+w]=(vec[x]-vec[x+1])/sqrt2;
        x+=2;
      }
      memcpy(vec,vecp,sizeof(T)*n);
      delete[] vecp;
    }
    static void Haar2Dhelper2(T *vec,int n,int span)
    {//one step vertical haar
      T *vecp=new T[n];
      const float sqrt2=sqrt(2.0);
      const int w=n/2;
      int x=0;
      for (int i=0; i<w; i++)
      {
        vecp[i] = (vec[x]+vec[x+span])/sqrt2;
        vecp[i+w] = (vec[x]-vec[x+span])/sqrt2;
        x += 2*span;
      }
      x=0;
      for (int i=0; i<n; i++)
      {
        vec[x]=vecp[i];
        x+=span;
      }
      delete[] vecp;
    }
    static void AntiHaar1D(T *vec, int n)
    {
      int w=1;
      T *vecp = new T[n];
      memset(vecp,0,sizeof(int)*n);
      const float sqrt2over2=sqrt(2.0)/2.0f;
      while(w<n)
      {
        for(int i=0; i<w; i++)
        {
          vecp[2*i]=(vec[i]+vec[i+w])*sqrt2over2;
          vecp[2*i+1]=(vec[i]-vec[i+w])*sqrt2over2;
        }
        memcpy(vec,vecp,sizeof(T)*w*2);
        w*=2;
      }
      delete[] vecp;
    }
    static void AntiHaar2D(T *vec,int n)
    {//Anti 2d haar on a nxn map
      int height=2;
      while(height<=n)
      {
        for (int i=0; i<height; i++)
          AntiHaar2Dhelper2(vec+i,height,n);
        for (int i=0; i<height; i++)
          AntiHaar2Dhelper1(vec+n*i,height);
        height*=2;
      }
    }
    static void AntiHaar2Dhelper1(T *vec,int n)
    {//one step horizontal anti-haar
      T *vecp=new T[n];
      const float sqrt2over2=sqrt(2.0)/2.0f;
      const int w=n/2;
      int x=0;
      for (int i=0; i<w; i++)
      {
        vecp[x]=(vec[i]+vec[i+w])*sqrt2over2;
        vecp[x+1]=(vec[i]-vec[i+w])*sqrt2over2;
        x+=2;
      }
      memcpy(vec,vecp,sizeof(T)*n);
      delete[] vecp;
    }
    static void AntiHaar2Dhelper2(T *vec,int n,int span)
    {//one step vertical anti-haar
      T *vecp=new T[n];
      const float sqrt2over2=sqrt(2.0)/2.0f;
      const int w=span*n/2;
      int x=0;
      for (int i=0; i<n; i+=2)
      {
        vecp[i]=(vec[x]+vec[x+w])*sqrt2over2;
        vecp[i+1]=(vec[x]-vec[x+w])*sqrt2over2;
        x+=span;
      }
      x=0;
      for (int i=0; i<n; i++)
      {
        vec[x]=vecp[i];
        x+=span;
      }
      delete[] vecp;
    }
  };

  typedef Cubemap<float> Cubemapf;
  typedef Cubemap<Vector3f> Cubemap3f;

  template<>
  bool Cubemap<Vector3f>::LoadFromPfm(const string &filename)
  {
    FILE *fp;
    fp=fopen(filename.c_str(), "r");
    if (fp==NULL)
      return false;

    char s[1000];
    s[0]=0;
    int width,height;
    int headerlen=0;
    if (fscanf(fp,"%[^\x0a]\x0a",s)!=1)
    {
      fclose(fp);
      return false;
    }
    if ((strlen(s)<2)||(s[0]!='P')||(s[1]!='f'&&s[1]!='F'))
    {
      fclose(fp);
      return false;
    }
    headerlen+=(int)strlen(s)+1;

    if (fscanf(fp,"%[^\x0a]\x0a",s)!=1)
    {
      fclose(fp);
      return false;
    }
    if (sscanf(s,"%d %d",&width,&height)!=2)
    {
      fclose(fp);
      return false;
    }
    if (width%3!=0||height%4!=0||(((width/3)&(width/3-1))!=0)||(((height/4)&(height/4-1))!=0)||width/3!=height/4)
    {
      fclose(fp);
      return false;
    }
    headerlen+=(int)strlen(s)+1;

    if (fscanf(fp,"%[^\x0a]\x0a",s)!=1)
    {
      fclose(fp);
      return false;
    }
    headerlen+=(int)strlen(s)+1;

    fclose(fp);
    fp=fopen(filename.c_str(), "rb");
    float *data;
    data=new float[width*height*3];

    if (fseek(fp,headerlen,SEEK_SET)!=0) return false;
    int ii=(int)fread((void *)data,1,width*height*3*(int)sizeof(float),fp);
    if (ii!=width*height*3*(int)sizeof(float))
    {
      delete[] data;
      data=0;
      fclose(fp);
      return false;
    }
    fclose(fp);

    cubesize=width/3;

    ChangeSize(cubesize);

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(POSX,j,i,cubesize);
      v[c.ToIndex()]=Vector3f(data+((cubesize*3-1-i)*width+cubesize*3-1-j)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(NEGX,j,i,cubesize);
      v[c.ToIndex()]=Vector3f(data+((cubesize*3-1-i)*width+cubesize-1-j)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(POSY,j,i,cubesize);
      v[c.ToIndex()]=Vector3f(data+((cubesize*3+i)*width+cubesize+j)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(NEGY,j,i,cubesize);
      v[c.ToIndex()]=Vector3f(data+((cubesize*1+i)*width+cubesize+j)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(POSZ,j,i,cubesize);
      v[c.ToIndex()]=Vector3f(data+((cubesize*0+i)*width+cubesize+j)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(NEGZ,j,i,cubesize);
      v[c.ToIndex()]=Vector3f(data+((cubesize*3-1-i)*width+cubesize*2-1-j)*3);
    }

    delete[] data;

    return true;
  }

  template<>
  bool Cubemap<Vector3f>::SaveToPfm(const string &filename)
  {
    char s1[1000];
    FILE *fp=fopen(filename.c_str(), "wb");
    if (fp==NULL) return false;
    sprintf(s1,"PF\x0a%d %d\x0a-1.000000\x0a",cubesize*3,cubesize*4);
    fwrite(s1,1,strlen(s1),fp);

    int width=cubesize*3;
    int height=cubesize*4;
    float *data;
    data=new float[width*height*3];
    memset(data,0,sizeof(float)*width*height*3);

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(POSX,j,i,cubesize);
      memcpy(data+((cubesize*3-1-i)*width+cubesize*3-1-j)*3,v[c.ToIndex()].Data(),sizeof(float)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(NEGX,j,i,cubesize);
      memcpy(data+((cubesize*3-1-i)*width+cubesize-1-j)*3,v[c.ToIndex()].Data(),sizeof(float)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(POSY,j,i,cubesize);
      memcpy(data+((cubesize*3+i)*width+cubesize+j)*3,v[c.ToIndex()].Data(),sizeof(float)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(NEGY,j,i,cubesize);
      memcpy(data+((cubesize*1+i)*width+cubesize+j)*3,v[c.ToIndex()].Data(),sizeof(float)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(POSZ,j,i,cubesize);
      memcpy(data+((cubesize*0+i)*width+cubesize+j)*3,v[c.ToIndex()].Data(),sizeof(float)*3);
    }

    for (zuint i=0; i<cubesize; i++) for (zuint j=0; j<cubesize; j++)
    {
      CubeCoord c(NEGZ,j,i,cubesize);
      memcpy(data+((cubesize*3-1-i)*width+cubesize*2-1-j)*3,v[c.ToIndex()].Data(),sizeof(float)*3);
    }

    fwrite(data,1,width*height*3*sizeof(float),fp);
    fclose(fp);

    delete[] data;
    return true;
  }
}
